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ABSTRACT 


The extended method of equivalent inclusion developed earlier by 
Fu (August, 1981}, (February, 1982) is applied to study the specific wave 
problems (i) the transmission of elastic waves in an infinite medium con- 
taining a layer of inhomogeneity, and (ii) the scattering of elastic waves 
in an infinite medium containing a perfect spherical inhoinogeneity. The 
eigenstrains are expai\ded as a geometric series and method of integration 
for the inhomogeneous Helmholtz operator given by Fu and Mura (1982) is 
adopted. 

The purpose of the study is to compare results, obtained by using 
limited number of terms in the sigenstrain expansion, with exact solutions 
for the layer problem and that for a perfect sphere. Two parameters are 
singled out for this comparison; the ratio of elastic moduli, f«(,\'+2y')/ 
(l+2y), and the ratio of the mass densities, h=p'/P* General trend for 
three different situations are showii: (i) A p»0, ( A A+2Ay)j^0, (ii) A pj^O, 

(AA+2Au)aO and (iii) A ,o?*0, (dX't2A y)^0. In the last case, material systems 
given in Truell, Elbaum and Chick (1969) are studied. 
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CHAPTER 1 


INTRODUCTION 

The scattering of elastic wav-f^s by an inhomogeneity embedded in an 
infinite homogeneous isotropic elastic medium has been studied by 
numerous investigators. An '’inhomogeneity" is a region in which dif- 
ferent material properties from its surrounding medium exist. Eshelby 
[1,2] developed the method of equivalent inclusion to determine the 
ji,elastic field of an ellipsoidal inclusion. An ’’inclusion” is considered 
to be a region which has the sajije geometric shape and dimension as the 
inhomogeneity but the same material properties as its surrounding medium 
after an eigenstrain is imposed within that region, Mai and Knopoff 
[3] appeared to be the first in applying Eshelby’ s result to form the 
scattering theory of a single sphere. Subsequently, Gubernatis [4] 
also used Eshelby *s result to study the long-wave scattering of elastic 
waves for an ellipsoidal inhomogeneity. In the above studies Eshelby 's 
solution for the static displacements was used as the first approxima- 
tion in the iteration. 

IVheeler and Mura [5] first applied the concept of eigenstrain to 
study composite materials, in which they considered the difference in 
elastic moduli between the inhomogeneity and matrix. Subsequently, 

Fu [6] presented a formulation for the elastodynaraics field of two 
ellipsoidal inhomogeneities embedded in an infinite elastic medium sub- 
jected to plane time-harmonic waves. He [7] later gave a complete 
formulation in extending the method of equivalent inclusion to dynamic 
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elasticity and gave some results for three-layered and five-layered 
media subjected to plane time -harmonic longitudinal waves. The scat- 
tering of plane waves by an ellipsoidal inhomogeneity is presented in 
[8.9]. 

The scattering of a plane compressional wave by a spherical inhotno- 
geneity in an infinite elastic medium has been studied by Ying and 
Truell [10] and by Pao and Mow [11] . They used the method of separation 
of variables to solve the wave equation which describes the incident, 
reflected (scattered) and refracted waves inside and outside the spheri- 
cal inhomogeneity. Because the solutions are expressed by a spherical 
coordinate and the inside and outside scattering field can match exactly 
along the boundary of the spherical inhomogeneity, the results by this 
method are considered as exact solutions. Some numerical results of 
the elastic scattering cross section were shown by Johnson and Truell 
[12], and some dynamical stress concentration factors around a spheri- 
cal cavity were found by Pao and Mow [14] . 

This research is concerned with numerical calculations according 
to the extended method of equivalent inclusion by Fu [6-9], and a com- 
parison with the exact solutions. Two cases are studied here. The 
first case is concerned with the three-layered problem and the other 
case is concerned with the spherical inhomogeneity problem. Results 
by the direct integration method for the three- layered problem and the 
method of separation of variables for the spherical inhomogeneity 
problem are compared with those obtained by the extended method of 
equivalent inclusion applied on the two cases. 


CHAPTER 2 


WAVE MOTION IN AN ELASTIC MEDIUM 


2.1, G overning Equation 

Consider the infinitesimal element of an elastic body with mass 
density p in the absence of body forces, the equations of motion are 

“jk,k ■ p 

where a dot indicates a differentiation with respect to time while a 
subscript comma indicates spatial differentiation. 

If the elastic material is linear and homogeneous, the stress- 
strain relation is 

a., » C., e C2) 

jk jkrs rs ^ ' 

where are the elastic constant. 

If the elastic material is isotropic, then the independent elastic 
constants are reduced to two and 


'jkrs 


= X 6 


jk rs 


y (5 . 6 , 

^ jr ks 


+ 6 


JS 


'kp 


(3) 


where X,y are Lame's constants and axe kronecker's delta, 

Substituting eqs. C2. 1.2) (2.1.3) into eqs. (2.1.1) the equations 
of motion in terms of the displacements are obtained as follows: 


C.uy) U 


5 


+ U 


U. 

1 




= P u, 


(4) 


or, in vector nutation as; 

(X+y) 77 • lT + yV^ = p if 
3 


( 5 ) 


4 


where 7 is the vector differential operator. 

The displacement vector can be decomposed as 

Ct’a7t|; + 7x5, 7*|a0. (6) 

where the '4> and ^ are the scalar and vector displacement potentials, 
respectively. 

A substitution of eqs. (2.1.6) into eqs. (2.1.5) leads to 
7[(X+2v)7^t|/ - p ii,] + 7 xfu 7^7 - p = 0 (7) 


and the !{i and 7 satisfy the wave equations 


2 1 " 


( 8 ) 


7^7 7 


(9) 


2 2 

;l.n which v^ = (X+2p)/p and v.^ = u/p are the velocities of longitudinal 
waves and shear waves, respectively. 


2.2. Initial and Boundary Conditions 

In the above section the general wave equations in the absence of 
body forces have been mentioned without considering external forces. 
Usually the external force can be treated as an internal or surface 
source (incident wave), or in the form of initial conditions (impulses) 
and boundary conditions (displacement or force or mixed boundary condi- 
tions) . Here the external force from sources, i.e. incident waves, will 
be discussed, 

If the incident waves are assumed to be plane waves, in general 
they are composed of compressional and shear waves. IVhen the incident 



5 


waves move in an elastic medium, the compressional and shear waves will 
propagate independently. But they can not travel independently if there 
is an inhomogeneity in the elastic medium. When the incident waves im- 
pinge on the inhomogeneity, compressional and shear waves will be re- 
flected back into the matrix while the same types of waves will be re- 
fracted into the inhomogeneity. Both the reflected and refracted waves 
must satisfy the general wave equations. If the elastic inhoitiogeneity 
is bounded to the matrix at all times, then the tractions and displace- 
ments must be continuous at the interface between the inhomogeneity and 
matrix. 


2.5. Time-harmonic Wave Problem 

If the incident waves are time -harmonic, then the reflected and re- 
fracted waves will also be time-harmonic waves of the same angular fre- 
quency, Therefore, the displacement of all the waves can be represented 
by 

0(r,t3 = U(r) exp(-io)t)e^ (1) 

where oj is the angular frequency and e^ is the unit vector in the direc- 
tion of wave propagation. 

The scalar and vector displacement potentials are 

i|;Cr, t) = ip(r) expC-iwt] (2) 

J(r,t3 - ({)(?) exp(-iwt)e^ C^) 

where e is the unit vector in the direction of shear wave. 

4 > 

Substituting eqs. (2.5.2) (2.3.3) into eqs . (2. 1.8) (2.1.9) , the 
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wave equations are roduoed to 

(7“ + d“) * 0 (4) 

(V‘ + e“) ^C?) « 0 (S) 

where a»w/Vj^, i5»w/v,j, are the wavenumber of longitudinal and shear waves^ 
respectively. 



CHAPTER 3 

INHOMOGENEITY PROBLEMS 


3.1, Transmission and reflection of waves in an infinite three-layered 
medium " ' ' ' ' ' ' 

A plane compressional incident wave which is simple harmonic is 
assumed to propagate in a three-layered medium, The geometry and ma- 
terial properties of the medium are shown in Fig. 1, The displacement 
of the incident wave is 

" '^o “ iwt)e^ (1) 

where is the amplitude and k, m are wavenumber and angular frequency ^ 
respectively, is the unit vector in the positive direction of z-a.xiSj 
i.e, plane wave propagation, 

3.1.1. Integration of the governing differential equation 

When the plane compressional incident wave impinges on the inter- 
face between the elastic inhomogeneity (medium II) and its surrounding 
medium (medium I), a compressional wave is reflected back into the 
medium 1, while a compressional wave is transmitted into region III [13]. 
Fox convenience of the ensuing discussion, the displacements and stresses 
associated with the incident, reflected, refracted and transmitted waves 
will be designated by the superscript (i), (r), (f) and (t) . Details 
for incident compressional plane waves along the +z axis are given below. 
The displacement of these waves are 

= U exp(ia - 2 - ioit) (1) 

Z O 'K 

« -Aexp(-iaTZ - iwt) 
z *■ y i 


( 2 ) 


8 


» B expCiajL^ - 

(33 

Cf3 „ u-^-^^Cz) expC-iwt) 

(43 


where A, B are the EunpUtude of the incident, reflected and trEins- 
mitted waves, respectively, , where the sub- 

script 1 denotes the matrix. 

The displacement function should satisfy the wave equation 


(V“ + ** ^ 


CS) 


where a 2 **^ ® " (x subscript 2 refers to the inhsimo- 

geneity , 


(f) 


For the simple case chosen, ^(z) are obtained as 


= Ccosa2Z + Dsina^z 


( 6 ) 


where C, D are constants. 

The stresses associated with the above displacements are 

= iaj^(Xj+2u^) expCia^^ - iwt) CO 

^22^ - -ia^CXj+2u^) A expC-ia^t - CS) 

= iai(Xi+2pj_3 B expCict^s - iwt) (9) 

ofP = a-CX,,+2vu)C-Gsina^z + Dcosao^) 

Tlie stresses and displacements must be continuous at the inter- 
face. Thus at 2 = - 5/2 the continuity conditions are 


1 

2 

,(r) , 

2 2 

cii: 

22 

„(x3 _ _(f3 

Of,„ - 0^^ 

CIS) 


OF POOR QUAtrrf, 
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and at z = 5/2 the continuity conditions are 


U 


(t) , uCf) 
z z 

,ct) , ^Cf) 


zz 


zz 


(13) 

(14) 


These continuity conditions give four simultaneous equations in 


terms of four unknown coefficients A, B, C and D. After solving the 
simultaneous equations, the unknowns A, B, C, D are found to be: 


cos- 


“2^ 


A=U exp(-ia, 6)[1 — r 


sin-^ 


'*2 

cos-j i/fh sin 


ct25 old a-5 

- sin— ^ +i/5h" cos- 


-] (15) 


cos 


“2^ 


B=U^expC-i=>i'S)[ ^ 


. “2« 
si n" ' 


2« ^2'^ ^2^ 
cos-^ ii/fh sin-y- sin— y cos-y 


•] 


(16) 


ot_5 ct-6 

C=U^expC-ia^6/2) / (cos— ^ i/fh sin-y) 


a^6 


a, 5 


D=~U^expC-iaj^5/2) / (sin— ~— + i/fF cos-y) 


(17) 

(18) 


where f = (X 2 + 2 y 2 )/ (^j^+2pj^) and h = P 2 /P 2 


The displacements and stresses in the whole medium can be obtained 
by using equations (1-4) and equations (^-10), respectively. 


3.1.2. Extended Method of Equivalent Inclusion 

The equations of motion for the inhomogeneity and the matrix are: 


"jk.k = »1 “j 

matrix 

(13 

*• 

'^jk,k ~ ^2 

inhomogeneity 

(2) 
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The stress-strain relations are; 


“jk = "rs 

'^jk " ^jkrs Ss 


matrix 


inhomo geneity 


(3) 

(4) 


where p^, and P2» are mass density and the elastic constants 

for the matrix and the inhomogeneity, respectively. 

Replacing the inhomogeneity with an equivalent inclusion which has 
the same material properties as the surrounding matrix after an eigen- 
strain is imposed in the inclusion, the governing equaticns are equations 
(1,3) and 

• fl 

inclusion (5) 


“jk,k - h 


“jk = S-tos 

= e - e* 
rs rs rs 


inclusion 


( 6 ) 

(7) 


where e , e® , e are the total strain, elastic strain and eigen- 

!]r 3 IT S 3T 3 

strain, respectively. 

The equivalence conditions derived and given in [7,9] are: 


= - AC,,. , inside the inclusion 


'jkrs ■'rs 


'jkrs rs 


'jkrs rs 


Ap iiW . = . Ap 

j jkrs rs,k j 


inside the inclusion 


where AC . Ap = Pj - Pj and 

cC13 / = c.'l’ 

jkrs ^rs,k jkrs rs,k jkrs rs,k 

U. = 

1 J 1 


(S) 

(9) 


( 10 ) 

( 11 ) 
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where the superscripts "m" and *'a" denote field that is induced by the 
presence of mis ‘-match and that is applied. 

Equations (3.1.2,5-3.1.2,11) are the general equations for the 

t 

equivalent inclusion with no restrictions on the geometric shape of the 
inclusion, i.e, these equations can be applied to an inclusion of arbi- 
trary geometry. 

For the three-layered problem, the equivalence conditions are re- 
duced to 

n=o 

^Apm- I 

n=o 

l [^i'"ho]A„ - - C^i - = -Ap«\ 

n=o 

and 

^(A,U2Ap) I = -CAU2Ap)E^ 

n=o 

^CAA-2Au3 I u'^h°lA„+'!'®(01BJ - CAi-2Pj)Ej = -(AA+2Ap3E^ (13) 


^(AX-2Ay) f t (.yEViDB^ = -(A.U2Ap)E^ 


2v1j^)Aq = -Apw 

+ 2yj^)A^ - -Apto^Hj^ (12) 
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where: == 2lo^ Uz' 

i^\ — X2 ~ 

Ay = p2 - Pj^ 

Cid,)" 

»m = — !- ‘>0 

E = — -i- U 

m ml o 

The eqs. (3.1.2.12) and (3.1.2,13) give 2(m+X) simultaneous equa- 
tionSi having matrix dimensions o£ 2(m.+l) x 2(m+l} if matrix represen- 
tation is used. 

The displacement and strain fields inside and outside the equiva- 
lent inclusion are; 


U^CZjO = U^(z) exp(-iwt) 


m • 

U,CZ) = I [A„*„Cz) * B 
n=o 


exp(-iwt) 


m i rt. 7 

= I * idiU^o“i 

n=o 


and the stress fields are 


a^jCz.t) - CT22C2) exp(-iwt) 


= ('^x + 2y^_) S2z^^^ 

^22(2) = + 2 y 5 ) e^^(z) 


matrix 


(14) 

(15) 

(16) 
(17) 


(18) 

(19) 

( 20 ) 


inclusion 
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Evaluation and integration of the 41 ^-integrals are shown in 
Appendix I. 

3.2. Scattering of waves by a spherical inhomogeneity 

A plane compress ional incident wave which is simple harmonic is 
assumed to move in the positive z-axis in an infinite elastic medium 
where an elastic spherical inhomogeneity is embedded. The configuration 
of the whole domain is shown in Fig. 2, The displacement of the inci- 
dent wave is 

- II exp(iaz - iwt)e ( 1 ) 

where is the amplitude and a, u> are the compressional wavenumber 

» 

and angular frequency^ respectively, and e^ is the unit vector in the 
positive direction of the z-axis. 

In this section two methods, the separation of variables ap- 
proaches and the extended equivalent inclusion method, are used to 
determine the displacements and stresses of the reflected (scattered) 
waves. Then, two measurable physical quantities, differential scat- 
tering cross section and the total scattering cross section far from 
the inhomogeneity, will be expressed in terms of the scattered asymp- 
totic values. The differential scattering cross section dP(w)/dQ [4] 
is a measure of the fraction of incident power scattered into a parti- 
cular direction, where dH is the differential element of solid angle. 

The total scattering cross section P (w) is the ratio of the average 
power flux scattered into all directions to the average intensity of 
the incident fields. 
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3.2.1. Method of separation of variables 

In this section, the work by Ying and ruell [10] and that by 
Pao and Mow [11] will be briefly introduced. 

When a plane compressional incident wave impinges on the surface 
of the elastic inhomogeneity, scattering occurs. Both compressional 
and shear waves are reflected back into the matrix (Medium I) while 
the same types of waves are refracted into the inhomogeneity (Medium II) . 
The potentials, displacements and stresses associated with the incident, 
reflected and refracted waves are denoted by the superscript (i), (r) 
and (f ) . The wave equations in terms of potentials are 


(V^ + = 0 

(V^ + = 0 

(V^ + a^^) = 0 

( 7 ^ + 0 ^^) = 0 


(7^ + = 0 

(7^ + B^”) = 0 

2 2 2 2 
Where /(Xj^+2y^) , /y^ 


02 = P 2 “ /^2 


= P2^ /(X2+2y2), 


( 1 ) 


( 2 ) 


(3) 


The displacement of the incident wave is 

exp(ia ^2 - iuit) 


(4) 


L 
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and the associated potential functions are 


'P 


Ci) 


5^ expCiajZ) 



» 0 


C5) 

( 6 ) 


can further be represented in terms of a spherical coordinate 
function 

't' = o' I CSn+m" jnCciir) P^(cos6) (7) 

1 n=o 


where 3-^s the spherical Bessel functions of the first kind, and 

P^(x) are Legendre polynomials. 

Because of axisymmetry of this problem, displacements, stresses 
and potential functions are independent of the spherical coordinate $ . 
Therefore, after solving the wave equations, the potentials of the re- 
flected wave are [12] ; 

' 00 

= I C-i)"^^ a(2n+l)Aj^ h^(a^r3Pj^Ccose) (8) 

n=o 


00 

tW , I a|-2n+l)B^ h_^(ejr)P^C=os6) CS) 

nso 

where the spherical Hankel functions of the first kind. 

For the refracted wave inside the inhomogeneity the potentials are 

[12] t 

CO 

= I a( 2 n+l)C^jj^Ca 2 r)P^(cos 0 ) ( 10 ) 

n=o 
n=o 


aC2n+i)D^j^(g2r3P^Ccos9) 


( 11 ) 
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The resultant waves in the medium I and medium II are 


l{l^ 5s 

(12) 

(j)^ = 

(13) 

t|^2 - 

(14) 


(15) 


The general representation of displacements and stresses in terms 
of the potential functions can be found in [15] : 




(16) 


II - 1. 

“ r 38 " 38 

rr 2 ^ r 3r 2 0^ r ^e''3r r-^ 

r 


_ _ 0,1 r B ~2ct X 1 . 1 n ,1, — >- 

^80 " 2 ^ r 3r 2 ^6'^ “ 2 


cote ^ 
36 


(17) 

(18) 
(19) 


cote 3 .n, .. In S-J), 
^ “ 30 - r ^03?^ 


, . s^-2a^ , 1 , cote 3* , °e'*' cote 3 „ „ 

7 T? * ~r M * ~2 - T~ 3e VJ 

r r 


( 20 ) 

( 21 ) 


13 13 3 

where = tr ~ and Dn = -rrror -sor(sin6 . For the incident field, 


r r 3r 


0 sin0 36 


'P - # 4>=0> ^”^1’ reflected field, tji = 4>=(p^^\ 

a=a^, 0=6j^ . For the refracted field, ct=ttp, 0=0.2 . 
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Consider the boundary conditions at r»a 
jjCi) + y(^) . y(f) 




= U 


r 

(f) 

e 


,Ci5 , , ,(f) 


rr 

,Ci) 


rr 

(r) 


o'Z' + a „ 

re re 


= 0 


rr 

(f) 

re 


(22) 

(23) 

(24) 

(25) 


Substituting eqs, (7-11) and eqs. (16-21) into the above continuity 
equations will form four simultaneous equations. After solving these 
equations j the unknowns will be obtained. Thus, the 

elastic field inside and outside the inhomogeneity can be determined. 

The total scattering cross section far from the inhomogeneity is 


[ 12 ]: 


p(u) = 4ira^ I (2n+l)[lA^|^ + n(n+l) |b^|^] 


(26) 


n=o 


and P(w) can be normalized as 


PW = 4 I (2n+l)[|A_^|^ + n(n+l) ^ 


(27) 


n=o 


3.2.2. Extended Method of Equivalent Inclusion 

The eqs. (3. 1.2.1 - 3.1.2.11) are general equations by this method 
and can be applied to the inhomogeneity problem of arbitrary geometry 
[7-9] . The equivalence conditions are recorded here as follows : 


+Ag = ” Apw Hg 


(1) 
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+A “ -iiPW 
sp sp 


( 2 ) 


+ ■57 2Au(djjj[0]A. + dstjkt“l''jk*'"*°stjk^'’^®jk*°stjkit°^®jki*-”^ 

* C4«stW2'‘l®st> “ ■C‘^“st=n» * 


. ^2Ap(d^,._p[0]A^*d3p.,_p[0]A.,.. . •+Dstjk,pt°l®jk*'’stjki.p['>5®Jki*- ■ •] 


where; Ap = P2 - Pj_ » " ^2 " ^ = VI2 " ^1 

H = U„ s=3 

s 0 

= 0 otherwise 

i 

H = iaU^ S“3 p=3 

sp o 

= 0 otherwise 

E = iaU„ m=3 n=3 

mn o 

= 0 otherwise 

E = -U m=3 n-3 p=3 

mnp o 

= 0 otherwise 



owgimw- 
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47IPjA.j[?] = + tjs - •i'js] 

* ♦kjs • ♦k.js’ 

- ^kit.skj * 

4ito to^d frl — ” *5. 'S, + '|). ■*‘1* •_ ,_j 

4TTPiU a . [rj jm ^,m :)n 


4iTpiW^dnjnjkt^] = - ^i + '^’k^m'^jn ^ ’^kjmn " '^kjmn^ 

‘*"f’l“ '^mnjk[r] = -VjBj [')' 

’ ^l“l^*,nm^jk 

‘*'''’l“^°mnijk'^’^^ ° ^*‘l^'^k,ijnm*'*’k,ijmn^“*'l®l f*k,in'’jm*‘*'k,iiil*3n^ 

where the <})- and ij;-integrals and their derivatives are evaluated for a 
sphere by using the method given in [8] and are listed in Appendix II. 
From Appendix II, £sj]^[0]» 

‘>5t3kit“]- -^stjk.ptO’ °stjk,pt°] 

(3. 2. 2.1 - 3. 2,2.4) can therefore be much simplified. 

If one assumes a two term expansion in the eigenstrains , then 

eqs. (3.2. 2.1) and eqs. (3. 2.2.4) will form 21 simultaneous equations 

while eqs. (3. 2, 2. 2) and eqs.. (3. 2. 2. 3) will form another 15 simultaneous 
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equations. The two sets of simultaneous equations are uncoupled. 

Solving the two set simultaneous equations for , B.. and 

1 3 K J K 

Bjkji’ displacements of the scattered wave can be written as [8,9] 


2 > '?> > > 
where: 4up,o) f. = -3i“<l)'S.+(b .-i 

•^1 m] 1 ^ mj ^,m] 


mj 


2 > 


2 > 


“’"=1“ ♦Aj * fn.rnj - 


mj 


- ^Pl^mkj * 




where; exp(ia R) 

♦ (?) = /// 5-^ =1V' 




= /// K 


exp Cictj^R) 


a 


R 


dV> 


r outside Q 


r outside Q 


expCiUj^R) 


'^k£ . . . s ^ m ^K* “^s R ^ outside Q, 


lkt?i = h;; 


r]^1 :z 

, kjl. . .s'- •' 3x^3x^. . . 3x^ 


tjj[r] 


exp(xBj^R) 

= /// dV 

fi 


R 


> ... expCiB R) 

9k [r] = /// X- K-^ dV’ 


R 


r outside 


f outside ^ 


H . , . s 
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[?] /// 


exp(i0^R) 


R 


dV* 


♦,kW 

S’ r*^T 3 > -•>, 

‘^,k)i. . “ 3Xj^9x^.V. x^ 


r outside 
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As |r| “► « the far field value of U can be expressed as [91: 

111 * " 


expia,r expi8-r 
+ D 

m m r m r 


C6) 


and the associated far field stress are 


expict,r 

^ C, a, 6 

nm 11 r k k mn 




expia^r 


fC a + C A 1 
^ m n n m’^ 


(7) 


+ 3. 


expi$^r 


(DA + D a n 
^ m n n m'^-' 


where a^ are direction cosines. 


The differential scattering cross section is therefore [9] 


dPCui) _ 
dn 


c . 

7 

a. 

D 

m 

X 1 

m 

U 

3, 

U 

0 

1 

0 


and the total scattering cross section is 


(8) 


dPCwj 


‘’'“5 = / dn 

and PCw) can be normalised as 


dQ 


C9) 


PM = dn 


ira 


CIO) 


Some important volume integral calculations are shown in Appendix II. 


CHAPTER 4 


COMPARISON OF COMPUTATIONAL RESULTS 


In this chapter numerical results are presented and compared for 
the three-layered and spherical inhomogeneity problems. Numerical 
results were made known by Truell and his co-workers for a perfect 
spherical inhomogeneity. In order to compare with these results, the 
same material properties are used and listed in Table 1. 


Table 1. Material Properties 


Material 

Compressional 
wave velocity 
(m/s) 

shear wave 
velocity 
(m/s) 

Mass density 
(g/cm3) 

Stainless Steel 

5790 

3100 

7.90 

Mg 

5770 

3050 

1.74 

A1 

6568 

3149 

2.70 

Ge 

5285 

3376 

5.36 

Polyethylene 

1950 

540 

0.90 

Be 

12890 

8880 

1.87 


^Vhen the extended equivalent inclusion method is used, there are 
complex matrices, with dimensions depending upon the dimensionless 
wavenumber and the differences in elastic constant and mass density 
between the inhomogeneity and the matrix, need to be solved. In some 
cases, the dimension of the matrix are so large that the numerical 
error can be relatively high. In order to get the best results, the 
IMSL (International Mathematical and Statistical Libraries) subroutine 
LEQ2C is used to solve matrices with double precision. The routine 
applies iterative improvement until the solution is accurate to machine 
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precision, If the matrix is too ill-conditioned to get effective itera- 
tive improvement, a terminal error is produced. 

4,1. Three -Layered Medium Problem 

The input parameters for both methods are the dimensionless wave- 
number 0^5, the relative ratio of elastic constants f and mass density 

h, where f is (\,+2ii.,)/ CX^+2u^) and h is 02 /^ 1 * displacement ampli- 
tude and the stress amplitude given in the preceding 

figures are nondimensionalized. 

The calculation of the exact solutions is simple and the dimension 
of the matrix is just 4x4. The calculation of the equivalent inclusion 
method is relatively complicated. The solution is an infinite series 
summation, and N^Cthe accepted number of terms to get convergent results) 
depends on f and h. In the calculation of the displacements and 
stresses from Eqs. (3.1.2.14 - 3.1.2.20), the sujnmation is considered 
to be acceptable until the ratio of the current term to the current 
partial term is less than 0.5%. 

In order to make detailed comparison, three cases are studied 
here. The first case considers the difference in elastic constants, 

i. e. f?4l and h=l. Fig. 3-7 and Fig. 8-12 display the displacement 

amplitude and stress amplitude vs respectively. The second case 
considers the difference in mass density, i.e. f=l and hj^l. Fig. 13-17 
and Fig. 18-22 display the displacement amplitude and stress amplitude 
vs a^(S. The last case considers the difference both in elastic 
constants and mass density. Fig. 23-27 and Fig. 28-32 display the 
displacement amplitude and stress amplitude vs From the results, 

it is found that the extended method of equivalent inclusion gives 
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excellent results which can be treated as the exact solutions. The 
values of in getting the convergent displacement and stress ampli- 
tude in the third case are listed in Table 2. 


Table 2. The Value of N£ for the Three-Layered Problem 



Ge in A1 

f=1.285 

h=1.985 

A1 in Ge 

f=0.778 

h=0.504 

Polyethylene in Be 
f=90.790 
h=2.078 

Mg in St 

f=0.219 

h=0.220 

St in Mg 

f=4.572 

h=4,540 

0.5 

6 

4 

6 

6 

6 

1.0 

6 

6 

6 

6 

8 

2.0 

8 

6 

6 

8 

8 

4.0 

12 

8 

8 

10 

12 

6.0 

16 

10 

8 

14 

14 

8.0 

18 

14 

8 

16 

18 

10.0 

22 

16 

10 

20 

20 


'4.2. Spherical Inhomogeneity Problem 

The input parameters are the dimensionless wave number the com- 

pressional and shear wave velocity and the mass density of the inhomo- 
geneity and the matrix. 

For the method of separation of variables, the results of the scat- 
tering cross section were plotted in Truell’s paper but the specific 

■j* 

values for different a^^a are not listed. 

For the equivalent inclusion method, there are two independent 
series summation. The first one is to get the inside function value 
in the eqs, (3.2.2. 1 - 3.2. 2.4), which are related to a^a and 
only. The second one is to get the outside function value in the eqs. 

(3. 2. 2. 5 - 3.2.2.10), which also depends on a^a and 

From the three-layered problem, it is known that the value of 

tThe data” sho;\rrT in this section are calculated from the computer program 
kindly supplied by Dr. J. Gubernatis. 
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depends on f and h. This is also true for the current problem, 
Besides, it is found that the value of N£ also depends on 
should be noted that the dimension of the matrix is proportional to the 
value of N£. Therefore, for large N£, the dimension of the matrix will 
be very large and the derivation and numerical computation to get the 
scattering cross section is very lengthy and time-consuming. Instead 
of finding the scattering cross section for large a^a, it is hoped to 
find the accepted value in what range of a^a when N£=l. Later, the 
result is called one-term solution when N£ is equal to 1. The eqs. 

(3. 2. 2.1 - 3. 2. 2.4) are reduced to: 

Ap«^[f^.[0]A. , , Aj = -Apoi-H^ (1) 

" 2Ap[d^^.[0]A. * 

* * 2^163^ =. , 2 AuE 3,3 (2) 


From the formulas in Appendix II, it is found that [0] 
d * . [0] are equal to zero, which make the eqs. (4.2. 2.1) and (4. 2, 2. 2) 

S tj 

uncoupled, i.e. the difference in mass density and the difference in 
elastic constants will have no coupled effects. After some manipula- 


tion, it is found that the nonzero variables are Aj, 3^2 and 

with B^^ equal to 3^2 symmetry. For low a^^a, the closed form 
solutions of these variables can be obtained and the closed form of the 


scattering cross section can also be obtained. 

Fig. 33-36 display the scattering cross section vs a^a from two 
method. It is found if a^a is less than 1, the tendency of the one 
term solution is good compared to the exact solution. For certain 
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material system the one-term solution represents » good approximation 
up to medium frequency range. The accuracy, however, is not good 
enough to produce relative errors less than 1% though goes down to 
0.01, The Fig, 37-40 display the comparison for the between 0,01 

and 0.1. It is therefore necessary to take more terms in the series 
for higher o^^a if this approach is to be used. 


CHAPTER 5 


DISCUSSION AND CONCLUSION 

It is worth studying the numerical characteristics of the ex- 
tended method of equivalent inclusion. By doing this, we will know 
what kind of inclusion this method can be applied to. Because of the 
simplicity and low computer time, the three-layered problem is chosen 
first for studying. Fig. 41-50 display the displacement and stress 
amplitude vs a^5, f and h, respectively. In these figures, the curves 
both for NCnumber of terms) = 5,12 and for the exact solution are 
shown. From these figures, it is found N£ will increase while 
or h increases or f decreases. Besides, the value of N£ for convergent 
stress amplitude is larger than the value for convergent displacement 
amplitude. It is also found that the displacement amplitude is less 
than the exact solution and the stress amplitude is larger than the 
exact solution at the larger a^5 if N is less than N£. 

Though the acceptable numerical results of the equivalent inclusion 
method in spherical inhomogeneity are not obtained in this study, the 
one-term solution supplies a very good tendency at low wavenumber. It 
is expected when the two-term solution (N=2) are applied, the accuracy 
at low dimensionless wavenumber will be very good and it will also 
supply good tendency in higher dimensionless wavenumber. And if N£ 
is large enough for different inclusion, the value from this method 
should be very close to the exact solution as in the three-layered 
problem. 
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From the above study, it is found will be large as the 
Cor or h or is very large or f is very small. Therefore, 

the dimension of the complex matrices will also be very large, which 
will probably cause the numerical singularity (algorthmic singularity) 
in high a^a range. The numerical singularity may be overcome after 
some numerical improvement skills are applied or other numerical methods 
are used. It is suggested that numerical scheme be carefully studied 
for large h, or small f at high wavenumber range. 
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EVALUATION OF SOME INTEGRALS FOR THE THREE-LAYERED PROBLEM 

I.l. Z** 

-h 


« C-D^'-^-Rj^Cah) 


for z > h 


” r=o (n-r)!(ia)^^'^ 




for -h _< z _< h 
for z < -h 


where: R^(ah) = J . 

r=o (n-r)!(ia)^ 




r=o (n-r) ! 


. o n +2 6/2 . ! ,1 

1.2. =^Z C2')V“I'’' Idl' 


- 6/2 


ia6 


•■= 1 0 e LnC«<5) 


y l^C-D^ nl Cctz)^ 1 

~ L O m J 


m=o 


(n-m) ! (iaz) 


m 


2 -^ 1.2 


n 

where : (a6 ) = ^ 


n 1 


.a 6 s 


n 


/ V I r ict 6 ^m 

m=o (n-m) ! ( ^) 


> , n-k+2 , (k) , , n-k+2 d . 

I.^. a 4 ^ " (z) = a — Y 

dz 


2 ^- "^-2 


- r C-l)^e^'^^C-l)"e 

^ ^ 2 

r\ — k “I , /* 1 t /* '\ n 

^ 14- (-1) n‘ (az) 

(n-m-k) ! (iaz) 


k_-iaz 


ia 6 


Ln(a 6 ) 


m=o 


n+k 


1 . 4 . [ 0 ] 
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2»0 


(1) both n and k axe not odd integer 

ia6 

^n-k+2^^(k)^Q^ » Lj^Co-5) 

icd 

- - L„C»S) - 2i] 

i 


(2) both n and k are odd integer 

ictS 

n-k+2, (k) _.k-2. 2 


a 


^^'.*^^[01 = - i'-'e Lj,C«5) 

jqo 

_ .k-2r^ 2 t 

= -1 [e Lj^CcuS) + — ] 

i 


(3) only one of n and k is odd integer 
^n-k+2^^(k) ^ Q 


.n_idZ 


A 1 \ ■t-VA *>* 

1.5. o"*>„Cj) = R.CciS) 

n 2i^ " 


_ e 


-laz 


2i2 " 


R„CiaS) 


1.6. R„Cc6) 


2i“ 

. C-il!®'"" 


n 


2i 


, R„Cc6 ) 


ia£ foS-^n 

, _ 2 P n! *•2 ^ 

where (cd ) = e I 


f ., , , ioS m 
m=o (n-niJ lC 
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for n < k 


for n > k 


for n < k 


for n > k 


for 2 ^ j 
for s < - ^ 


for z > 


for 2 < - 


m 1 o » io |<>' 
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APPENDIX II 

EVALUATION OF SOME VOLUME INTEGRALS FOR 
THE SPHERICAL INHOMOGENEITY PROBLEM 


From Ref. [6,8], we obtain; 
II. 1. V(r) s /// dV» 

n ■ 


r 


y y "y^ y 

Jo 1=0 k=o 


.% . R n-i-k 
z 


/// 


n 3x'^ 3y'*^ 3z' 


dx* dy* dz* 


B J x_x_ . . . X. 


n-o 


nl p q 


nth power of x 


0> 


[ I 

m=o 


/ « N in 2ni / \ 
C-l) g rCn),m 


C2m)l '-pq...v * ^ ”pq. 

where: R = jr - r' | 

;Cn),ni fjj — ? dV’ 

nn. .V Jjj Pv^ j . . .ayl ^ 


i I 




-pq...V ' 3X^3X^...ax; 

3Cn),m ^ III dV 

p(r') = 1 


II. 2. If/., 


Ic^ « • • s 


(?) 




B ^ 144_x_x_...x 

n=o 


n! "p^q' ■ '”v 

nth power of x 


00 . ,sm 2m . s 00 

_L C2m) I pq- • -V i, 

m=o ^ m=l 


y 2m- 

C“lj a 

C2m-1] ! 


inside 0 


,m , 

. .V •' 


inside 


i 1 

pq. . .v^ 


31 
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where ; 


-Cn),m 

pq . . . V 


gCn),m 

pq--*v 


* /// PC?') 

a /// p(r') 


3x^3X^..,3x; 


dV‘ 


dV 


p (r ' ) ® ‘ 
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In what follows we develop these formula specially for a sphere of 
radius a: 


II. 3. 


and s**^^"*^ for p(r) = 1 and fi is sphere 
pq...v pq...v 


CD 


cCo)»n» 

sCo),m 


= 4H 


a 4H 


„2m+2 

a 

2m+2 

2m+l 
a 

2m+l 


C2) 



S 


C2),m 

pq 


0 for m=0 


C3) 


cC4),m ^ 
pquv 


0 for m=0,l 
(2m-l) (2m-3)2m 

C2m-1) (2m-3)2m 


4n 
S ' 

4IT 

15 


p=q=su=v 

for any two equal pairs 
of indices 


otherwise 


;C4),m ^ 
pquv 


C2m-2)(2m-4)C2m-l3a^’”“‘^ ~ 
C2m-2) f2m-4) (2Tn-l) 


p=q=u=v 

for any two equal pairs 
of indices 


0 


otherwise 
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_ ^(3),ni _ _ _(3),in _ - 

p " pqu " p " pqu “ 


pqu . . . V 


and S 


(n),m 
'pq. . . V 


for p(r') = x| and is sphere 


p(l),m (2m-l) ^2m+2 4II . 
p “ 2m+2 ^ ■ 3 *^pi 

5;(l),ni „ 2m-2 ^2m+l 4H 
p “ 2m+l 3 pi 


rC3),m _ 

pqu ■ 


0 for m=0 

(2m-13 C2m-3)a^"^ 

(2m-13 (2m-3)a^‘‘‘ 


r- 

4n for any two equal pairs of 
15 indices 

0 otherwise 


^pqu’” = 

= (2m-2) C2m-4) 


2m- 1 

4K 

P=q=u=i 


5 

2m-l 

4n 

for any 


15 

indices 


two equal pairs of 


0 otherwise 


cC5)>m ^ 
pquvw 


0 for m=0, 1 

C2m-1) (2m-3) (2m-5)2m a^”*"^ ^ D 


D 

D 


4n 

i=p=q=u=v=w 


7 


4n 

/i=p q=u=v=w 

p,q,u,v,w can change the 

35 

tp=q i=u=v=w 

order arbitrarily 


D = 

D = 0 


4n 

105 


i=p q=u v=w 
otherwise 
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®pqivw ° C2m-2)(2m-4)C2m-6)C2ni-l)a 


2m-3 
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D = 

4n 


7 

H 

II 

II 

.Q 

II 

C 

II 

< 

II 

D = 

4H 

i=p q=u=V=w 


35 

1, p=q i=u=v=w 

D = 

4n 

105 

i=p q=u v=w 

D = 

0 

otherwise 


Co),m ^ 5 jCo),m ^ p(2),m ^ 5(2),m ^ (,(4),m _ g(4),m 

pquv ~ pquv 


(4) = S 

pq pq 


= 0 


II. 5. 




m-l 2m-l 


(o) ,m 


m=o 


mil 


f-’i * roi - f 11^ r"y (- 1 ]”°^°' I i y s(2).m, 

C2) ♦,pq[0] - C-13 *1^1^ (2^.1), ] 


(2m-l) I pq 
,m-l .2m-l 


r-i) * rni = r.n'*r y C-dV” p(4).m .. “ (-ir"-''c.^'"-% C4].m. 

^ ^ ^’pquv^ ] C ) (2m)! pquv (2m^l) ! pquv ^ 


^^pCO] = 'i^^pqutOJ = ^>pquvwf°^ = ° 


II. 6. 


kr) 


J J J 


r outside Q 




K 


“ n n-2. , ..n 

y T y llii 

^ Z\ k'Cn-Jl-k)! 

n=o ^=0 k=D 


atiC^^E^ 

^xhy^dz^-^-^ 


///Cx')^(y’3^(2’)”'^"^ dx'dy'dz' 


for |r| 

-> » 

and n is 

sphere 


Op 

n 

n-.1 


J(r) = 

I 

1 

z 

FCn,1,k,a,a,r) 


n=o . 1 = 0,2 k=o,2 


35 


where : 


II. 7. 


FCn.«,k.a.a,i) = C-io)" HElHill )’'(e^)" 


-A-k 


ITT ’ri+3 
4IIa 

(n+3) (n+1) 


Hi 

2 ' 


(n0c|0(|oc^^^0 


e ,e , e are direction' cosines 

A y i 


f'mW = 37- fW 
m 


= STsTT'I'W 

m j 


= I . .1 . I ’ P(n,^.,k,a,a,r) 


for |rj -> « and Q is sphere 
“ n n-i 

I I I 

n- 0,2 S-=o,2 k=o,2 

> ^ « n n-A 2 

= I I ^ (ia) e e. • F(n, A,k,a,a,r) 
■' n=o,2 ^-= 0,2 k=o, 2 •' 
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Fig. L Displacement amplitude as a function of,c»i6for the 
three-layered problem, A1 in Ge, with h*1.0 




ORIGINAL PAGE IS 
OF POOR QUALITY 

-EXfiCT SOLUTION 

©CURRENT SOLUTION RT Z/6 =-0.5 


o ^CURRENT SOLUTION RT Z/<5 = 0.5 



iFig. 5 Displacement amplitude as a function of 6 for the 
three-layered prohlem .Be in Polythylene, with h=*1.0 


ORIGINAL PAGE IS 
OF POOR QUALITY 

“EXACT SOLUTION 

oCURRENT SOLUTION AT Z/« ^-0.5 
-s^CURRENT SOLUTION AT Z/<5 = 0.5 



DIMENSIONLESS WAVENUMBER oj 6 


Pig. 6 Displacement amplitude as a function of a 5 for the 
three-layered problem. Mg in Stainless Steel, 
with hal.O 












ORIGINAL PAGF. IS 
OF POOR QUA^:TY 


-EXfiCT SOLUTION 

oCURRENT SOLUTION AT Z/5 =-0.5 


o J.CURREMT SOLUTION AT Z/i *0.5 

O 



Fig. 9 Stress amplitude as a function of 5 for the three- 
layered problem, A1 in Ge, vith h=1.0 





ORIGINAL PAGE IS 
OF POOR QUALITY 


-EXACT SOLUTION 

©CURRENT SOLUTION AT Z/ 6 =-0.5 


o .^CURRENT SOLUTION AT Z/ <5 = 0.5 



Fig. 11 Stress amplitude as a function of 6 for the three 
layered problem. Mg in Stainless Steel, with h=1.0 



ORIGINAL PAGE IS 
OF POOR QUALITY 

-EXACT SOLUTION 


^CURRENT SOLUTION AT Z/fi=~0.5 
a ^CURRENT SOLUTION AT Z/ 6 = 0.5 

a 



DIMENSIONLESS WAVENUMBER <5 


Fig. 12 Stress amplitude as a function of a 5 for the three 
layered uroblem , Stainless Steel in Mg, vi%h hsl.O 























ORIGINAL PAGE 13 
OF POOR QUALITY 

-EXACT SOLUTION 

©CURRENT SOLUTION AT Z/ 6=^-0. 5 
o acurrent solution at Z/ 5 = 0.5 

o 


1 



dimensionless WRVENUMBER ai5 


Fig. 20 Stress amplitude as a function of ai 6 for the three 
layered problem, Be in Polythylene, with f^l.O 


OniQlNAL PAGE IS 
OF POOR QUALITY 


-EXACT SOLUTION 

©CURRENT SOLUTION AT Z/<5=-0.5 


o ^CURRENT SOLUTION AT Z/ <5 = 0.5 



Fig. 21 Stress amplitude as a function of 5 for the three- 
layered problem. Mg in Stainless Steel, vith f»1.0 


ORIGINAL PAGE IS ' 
OF POOR QUALITY 


-EXACT SOLUTION 

©CURRENT SOLUTION AT Z/6--0.5 


o ^CURRENT SOLUTION AT Z/ fi « Q.5 



DIMENSIONLESS WAVENUMBER « 


Fig. 22 Stress amplitude as a function of for the three 

layered problem. Stainless Steel in Mg, with f=1.0 











DISPLfiCEf-fENT AMPLITUDE 

J.OO 0.10 0.20 0.30 O.iiO 0.50 0.60 0.70 O.SO 0.90 1.00 

I I I ! I I ♦ f I I ^ 


O^INAL PAGs IQ 
OF POOR QUArnv 


-EXACT SOLUTION 

©CURRENT SOLUTION AT Z/6 =-0.5 
^CURRENT SOLUTION AT Z/6 = 0,5 



Fig. 25 Displacement amplitude as a function of 6 
the three-layered problem, Be in Folythylene 










OWGINAt PAGE 13 
OF POOR QUALITY 


-EXRCT SOLUTION 

©CURRENT SOLUTION AT Z/6=-0.5 


a ^.CURRENT SOLUTION AT 2/ = 0.5 



DIMENSIONLESS WAVENUMBER 

Stress amplitude as a function of <5 for the 
three-layered problem, Ge in A1 

































ORIGINAL PAGE IS 
OF POOR QUALITY 


75 


o 



?igV 38 Scattering cross section as a function of lov 

for the spherical inhomogeneity problem ,A1 in Ge 






-EXfiCT SOLUTION 
©CURRENT SOLUTION M 
A CURRENT SOLUTION N 


; 

a 


CD 


o" 


OiiMENSIQNLESS WRVENUMBER qa 

Fig. Displacement amplitude at the point zjS *0.5 as 
a function of a j ^ for the three-layered problem 
Stainless Steel in Mg 







ORIGINAL PAGE IS 
OF POOR QUALITY 

-EXACT SOLUTION 
©CURRENT SOLUTION N =6 


o ^CURRENT SOLUTION N =12 
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Pig. 43 Displacement amplitude at the point z/<f»0*5 as a function 
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Fig. 41» Displacement amplitude at the point z/S’»0.5 as a 

function of f for the three-layered problem a the 
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Stress amplitude at the point z/^=0.5 as a function of f 
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Fig, Displacement ainplitude at the point z/S-0.'? as a function 
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Fig. 48 Displacement amplitude at the point z/<T =0.5 as a function 
of n for the three-layered problem at the wavenumber 
(/. =10.0 and f=1.0 
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Fig. L9 Stress amplitude at the point »0.5 as a function 
of h for the three-layered problem at the wavenumber 
^*2.0 and f »1 . 0 
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